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ABSTRACT 


The longitudinal stability of an aircraft at or near 
stall was examined using the digital computer as an experi- 
mental tool to solve the longitudinal equations of motion. 
A linear analysis determined the effect of lift curve slope 
variation. An investigation was made to identify the non- 
linear lift curve variations needed to create the often 
abseryed «"rocking=chair"™ ‘or | porpoising” stall- trait.. ~ The 


characteristics of this limit-cycle oscillation were examined. 


7 oN ee 


: ; a ‘ aa Nid 
4 sinha 
sain? i. in aaa. 
a4 noes - 1p eee a “ 


TABLE OF CONTENTS 


I. INTRODUCTION----------------------------------------- 9 
II. AIRCRAFT MODEL AND ASSUMPTICNS----------------------- 11 
III. DISCUSSION AND RESULTS------------------------------- bh 
IV. CONCLUSIONS------------------------------------------ 20 
APPENDIX A. TABLES------------------------------------- 22 
APPENDIX B. FIGURES------------------------------------ 23 
APPENDIX C. COMPUTER PROGRAM 2D OUTPUT---------------- 37 
APPENDIX D. MODIFICATIONS TO COMPUTER PROGRAM---------- 68 
LIST OF REFERENCES---------------------------------------- 70 
BIBLIOGRAPHY---------------------------------------------- TL 
INITIAL DISTRIBUTION LIST--------------------------------- 72 
FORM DD 1473---------------------------------------------- 73 


a 


: an en eae 


ws se 
oe ey . in eR: aan 


‘ 


; na si a 


: _— Bipasha en we ca inmost ee ane. TT 


= 2 : 


a - 
5; i aren Tia Ni 
a . Pom . | ° . 
: 


TABLE I. 


LIST OF TABLES 


CHARACTERISTIC EQUATION ROOTS 


FIGURE 


LIST OF FIGURES 


Aircraft Axis System and Coordinate Transformation--23 


Root Locus Plot for varying CK SS 24 


Damped Oscillation for Impulse Elevator Input 


a. Lift Coefficient versus Angle of Attack--------- 25 
b. Aircraft Velocity versus Time------------------- 26 
c. Aircraft Angle of Attack versus Time------------ 26 
d. Aircraft Lift Coefficient versus Time----------- 27 
e. Aircraft Pitch Angle versus Time---------------- 27 


Aarerart Lift Coefficient versus Angle of Attack 
with aC, (1G SS a SS 28 


Limit Cycle Oscillation for Impulse Elevator Input 


Bo ees US a SS 29 
b. Aircraft Lift Coefficient versus Time----------- 30 
c. Aircraft Angle of Attack versus Time------------ 30 
d. Aircraft Pitch Angle versus Time---------------- a 
e. Aircraft Velocity versus Time------------------- 31 
Pitch Angle versus Time for varying AC, ee 32 


Pitch Angle versus Time for varying Elevator 


Limit Cycle Oscillation for Step Elevator Input 
a. Aircraft Angle of Attack versus Time------------ 34 


b. Aircraft Pitch Angle versus Time---------------- 34 


Ab 7 ie ... ; ty Le os : @ TA5 rt, 
? ' 
<-> 7 : de ' we -* 
Mp: Vatkes Th 7h 
j >, ) nv _ ' ce : 
} 
7 7 we oa oe 
< > ia 5 - oe _ rn 4 » os 
‘ _ so 
: ps ae a _ o ® _ i Ys 7 - 
Cy 
— 


i ie Se ei 7 mee MA a «4 | 
a5 i. = es 7 eh pil Ree as Pree: ee ne ei. ) 
a 


dows = a peru 7 


2 Pere Ss Pye Ses wet Peabo ts banag 2 a 
etal sentephegam ees © +; 
ae —— ieieatnaas ers ae 


Peace x Peer Bec occa ne 


juqnT 10develt qe22 tot aOigelitoro sic tide Se < 
bE ------------on IT eversy doesté Yo alpnA 23etsaik 4a a 


$f -----------~-----gnit urexe¥ alyes dears dieertip .d 


10. 


Limit Cycle Oscillation for Ramp Elevator Input 
a. Aircraft Angle of Attack versus Time----------- 
b. Aircraft Pitch Angle versus Time--------------- 


Limit Cycle Oscillation for Stick Push-over 
from Aircraft Stalled Condition 


a. Aircraft Angle of Attack versus Time----------- 


b. Aircraft Pitch Angle versus Time--------------- 


Te 4 7 
a] oo fi 
My wm Wy) @ 
le * 
re 
“ y 
| 
- 
cD 
> 
ad ® f 
= be? “as. 
_ 
2 elo 
a he 
i 
7 7 Pie _ 
7 > "Th _ 7 
_ 7 7 : 4 : : 
2 . so 
=> 
_ an 


7 eo 
“ae By nen eee . 


TABLE OF SYMBOLS 


Drag Coefficient 
Lift Coefficient 


Lift Coefficient derivative for 
angle of attack per radian 


Maximum attainable Lift Coefficient 
Pitching Moment Coefficient 


Pitching Moment Coefficient derivative 
for angle of attack per radian 


Aerodynamic Force Coefficient derivative 
for angle of attack rate per radia per 
second 


Change in Lift Coefficient on the Lift 


Curve 
Time seconds 
Aircraft velocity in x direction feet per second 
Aircraft resultant velocity feet’ per second 
Aircraft velocity in z direction feet per second 
Angle of Attack radians 
Elevator deflection radians 


Root of the characteristic 


equation per second 
Undamped natural frequency radians per 
seconds 
Pitch Angle radians 
Pitch Angle rate of change with radians per 
time second 


Damping ratio 


2bao00e2 
Pas ifear 


Raq: ers ibex 
Diol ae 


Yéeq =@nhelbsat. 


fjiw sprsis io tex rips 


Olja: pebqest i | 


ACKNOWLEDGMENTS 


The author wishes to thank Professor Louis V. Schmidt, 
Professor of Aeronautics at the Naval Postgraduate School, 
for his guidance and support as thesis advisor. The author 
would also like to express his appreciation to Lieutenant 
Carlton W. Saul for his assistance in the preparation of the 


computer program used in this thesis. 


Lt. INTRODUCTION 


Federal and military aviation requirements stipulate that 
all aircraft possess stalling characteristics which comply 
with established criteria. The Federal Aviation Regulations, 
Part 23 of "Airworthiness Standards" state, in part, that 
"acceptable stalling characteristics be demonstrated...in 
straight flight with wings level...where the primary control 
manipulation is a steady progressive upward movement of the 
elevator until the aircraft is stalled." A desirable stall 
characteristic ina longitudinal sense would include a small 
but distinguishable drop in the nose of the aircraft, which 
can not be overridden by the pilot. Also, if the elevator 
control were eased forward, the aircraft should promptly 
return to unstalled flight. 

The exact definition of aircraft stall will differ among 
the pilot, wing design specialist, stability and control ex- 
pert, and all others who will have their own interpretation. 
Highly sweptwinged, supersonic aircraft have discredited the 


idea of flow separation, C and a stable nose-down 


L max' 
pitching moment all occuring simultaneously. The Federal 
Aviation Regulations, Part 25, defines the calibrated stalling 
speed for a civil air transport as"...the minimum steady speed, 
in ‘knots, sat whichsathe aircraft is.icontrolhableud:"»iFerathe 
purpose of this treatise, stall may be described as the con- 


dition when flow breakdown, primarily over the main wing, 


causes Significant nonlinear effects in the moment and lift 
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characteristics. The effect of an asymmetric flow condition 
leading to wing drop will not be considered, since only the 
longitudinal stability traits are being analyzed. 

A considerable amount of literature is available on the 
Subject of aircraft static stability at stall, and it de- 
scribes the aircraft's initial tendencies, or stall charac- 
teristics. Some literature is available concerning dynamic 
stability at deep stall, or large angles of attack. There 
is, however, very little treatment of airplane dynamics in 


the area of C of post-stall time histories in that 


L max' 
Marrow region of the stall break, or of the oscillations 
Ghat might occur, 

Most pilots with any degree of proficiency have ex- 
perienced a "rocking-chair" or "porpoising" type of stall 
where, with constant elevator control, the aircraft's nose 
oscillates vertically about some point as the wing stalls, 
recovers, rotates back up into the stall, and so on. It 
was the purpose of this study, therefore, to numerically 
investigate the ee eee required to create this type of 
limit-cycle, or possibly divergent, rocking-chair motion, 
using the complete set of aircraft equations of motion, and 
to attempt to describe the major factors influencing this 
phenomenon. | 

Because of time considerations and the availability of 
data, this study was limited to only the longitudinal sta- 
bility, with no roll, yaw, or coupling considered, of a 
small, straight-winged jet airefLart, with zero thrust. All 


computations were accomplished at the W.R. Church Computer 
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Center, Naval Postgraduate School, Monterey, California, 
using the IBM 360/67 digital system and standard FORTRAN 
language. The study was accomplished during the period 


September 1971 through June 1972. 


The input data used in this tae tnese Lor an F-94 
Barerart obtained from Ref. I. The aircraft is a straight- 
winged, tandem-seated, single-engine jet, with the center of 
gravity positioned well within the fore and aft limits. A 
zero thrust landing configuration at sea level on a standard 
day exists throughout this analysis. The aircraft has 
basically linear noble iad S characteristics except for lift, 
drag and pitching moment coefficients, which are all functions 
of angle of attack. Terms relating to such variables as 
elevator control effectiveness and pitch and angle of attack 
damping are also assumed to be linear. The aircraft was 
considered stick-fixed, with the elevator being the only means 
ef longitudinal’ control and* trim.” The case. of pitch-up at 
stall was neglected. Figure 1 shows the aircraft axis system 


and the related variables. 


III. DISCUSSION AND RESULTS 


Static longitudinal stability may be considered as the 
tendency of the aircraft to return to static equilibrium, 
while dynamic longitudinal stability considers the resulting 
aircraft motion as a function of time. The existence of 


Static stability does not necessarily intimate the existence 
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of dynamic Stability. However, the existence of asymptotic 
dynamic stability in the sense of Lyapunov (Ref. 4) implies 
Static stability. 

The factors influencing aircraft stability at a stalled 
condition must be known in order to predict or create an 
instability. The airplane has six degrees of freedom: ro- 
tation in roll, pitch and yaw, and translation in the 
horizontal, vertical and lateral directions. The principal 
variables in the longitudinal motion are: 

(1) the pitch attitude of the airplane, 

(2) the angle of attack, 

(3) the flight velocity. 
The equations of motion describing the longitudinal motion of 
the aircraft are found in Ref. 2. These equations are or- 
dinary nonlinear differential equations in five variables, 
with time as the independent variable. A complete derivation 
of the equations and their application to the problem of 
Bmarerart Stall may be found in Ref. 3. 

To perform a sensitivity analysis on the stability de- 
rivatives near stall, the equations were simplified by as- 
suming a constant elevator position and a horizontal initial 
flight path. Reference 2, Chapter 6, reduces the equations 
to three homogeneous algebraic equations in the unknowns 
a, ; and 6. By setting the determinant of the coeffi- 
cients equal to zero, the roots, A, of the characteristic 
equation may be found. Expansion of this “stability deter- 


minant" results in a quartic equation for A of the form: 
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The solution of this equation results in two pairs of complex 
conjugate roots. These es tig MS be expressed either in non- 
dimensional time or real time, with our interest focusing 
on E, the damping ratio, and Ow the natural frequency. 
These numbers are easily calculated along with the resulting 
period of oscillation. For a dynamically stable aircraft, 
the real part of all of the roots must be negative; that is, 
the damping ratio is positive. An instability results if 
any of the damping terms are negative. 

Upon analysis of the equations of motion, the only sta- 
bility derivatives making a significant contribution to 


Stability were: C ‘es , and Clk By varying these 


Ix 7 “me 
parameters and observing the damping ratios of the resulting 
roots of the stability determinant, a feel for the important 
factors in aircraft stability may be developed. 

Only two stability derivatives are of sufficient magnitude . 


£O Contribute to dynamic stability. The first of these is 


Cc 


y the static stability derivative. As C becomes 
mex Mox< 


positive, the aircraft becomes statically unstable. This 
condition will not be considered since most aircraft demon- 
Strate positive static stability through the stall, and dy- 
namic stability is of concern here. The second stability 


derivative of concern as ¢ the lift curve slope. The 


LO ’ 


results of varying lift curve slope are shown in Table l. 
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The aircraft demonstrates dynamic stability when rad is 
positive, zero or slightly negative ..“When the 1ift curve 
slope reaches approximately -5.0, one pair of the roots of 
the characteristic equation has a positive real part which 
indicates aircraft instability. The instability remains as 
the lift curve slope is increased negatively (c.f., Fig. 2), 
A negative slope of the lift curve could occur after the 
Maximum lift is reached in the area of stall, or in the case 
of a sharp stall break occurring locally. 

It must be noted at this point that only the short period 
motion was affected by this analysis. One would expect to 
find, therefore, that the aircraft has a stable long period 
or phugoid mode, but is unstable in angle of attack, pitch 
angle and normal acceleration (g's), with the velocity re- 
maining relatively constant as is characteristic of the short 
period mode. 

The computer program used to mathematically "fly" the 
aircraft is discussed at length in Ref. 3. Basically, there 
are four variables affecting the longitudinal motion: the 
resultant aircraft velocity, V, the pitch ang wey OQ, the rate 
of change of pitch angle, 6, and the angle of attack, &. 
All are functions of time, t. The rate of change of angle 
of attack is considered small and is neglected. Five dif- 
ferential equations were used to describe the longitudinal 
motion, the fifth occurring by resolving the resultant air- 
craft velocity, V, into its components, u and w, in the x and 
z directions, respectively. The equations were solved by a 


fourth-order Runge-Kutta integration scheme. Input data 
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consist of the aircraft parameters such as wing area, moments 
of inertia, thrust coefficient, air density, etc., anda 
table of data with aircraft lift, drag and pitching moment 
coefficients as functions of angle of attack. These input 
values remain constant throughout any maneuvering the air- 
craft is required to perform. The only initial conditions 
read into the program are the see wee gross weight, desired 
angle of attack and angle of pitch at trim. 

Subroutine. SPLIN1 calculates the value of Ch, ae and eS 
for the angle of attack specified from the data table by 
interpolation using a cubic curve-fitting scheme. Subroutine 
TRIM uses these data for a calculation of the corresponding 
velocities, elevator angle and thrust using the force and 
moment equations at static equilibrium. These calculated 
data are subsequently used as the initial conditions to start 
the Runge-Kutta integration scheme. Aircraft maneuvering 
and trim is performed by logic statements controlling elevator 
angle. Appendix C presents the entire computer program. 

The effect of a large negative value of lift curve slope 
(CL y = 10.0 per radian) beyond the stall was numerically 
introduced in the table look-up data for angles of attack 
greater than 0.4084 radians. Angles of attack escuchar this 
corresponded to the basic airplane situation as shown in 
Figure 3a. For angles of attack greater than this value, the 
lift curve slope was linear with a negative slope of -10.0 
pert fadian. “This type of drop-off in lift coefficient 1s° not 


unlike that encountered in wing leading edge flow separation 


OL Stall, 
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The aircraft was flown (numerically on the computer) at 
an initial trim condition of & =0.4083 radians (very near 
Ch a 9 =0.0989 radians, and a gross weight of 12,359 lbs. 
for three seconds. Then a trailing-edge-up elevator impulse 
(back stick) of 0.50 degrees (0.008725 radians) was inserted 
for 0.05 seconds. The resulting time histories show. the 
aircraft to exhibit divergent oscillations, as predicted 
in Fig 2, in angle of Sereek , pitch angle, and lift coef- 
ficient. The aircraft velocity also shows divergence, al- 
though the amplitude was very small compared with the other 
parameters, thus confirming the short period tendencies. 

The case of stall with a smooth lift curve as defined 
by Fig. 3a was investigated. Aircraft time histories with 
the same initial conditions are shown in Figures 3b through 
3e. This type of lift curve is typical of a wing trailing 
edge flow separation or stall. The maximum negative slope 
reached was ot =-2.0. The damping ratio for hae case 
was positive, with negative real parts for the roots of the 
characteristic equation, indicating an asymptotic type of 
dynamic stability in the sense of Lyapunov, with a dying out 
of the oscillatory motion to a static equilibrium condition. 
Mathematically, the aircraft model and the program used to 
Manipulate it have demonstrated both dynamic stability and 
instability in the region of stall based on the slope of the 
lift curve at angles of attack beyond that corresponding to 
Maximum lift. 

The problem of "rocking-chair" stall would indicate a 


limit-cycle type of oscillation. Reference 4 describes a 
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limit cycle as a nonlinear oscillation..."due to the presence 
of nonlinear terms in the differential equation." A stable 
limit cycle is characterized as being a unique closed curve 
in the phase plane, to which all other nonclosed trajectories 
approach in the form of spirals winding onto the limit cycle. 
It is a self-sustained oscillation, independent of the initial 
conditions. The differential equations governing the longi- 
tudinal motion are nonlinear as well as the lift, drag and 
moment coefficients being functions of angle of attack. With 
an instability in the system at static equilibrium due to the 
negative lift curve slope, it is possible for a limit cycle 
to exist. 

To prevent an oscillatory divergence as shown previously, 
the instability was bounded by a stable lift curve slope. It 
may be likened to a reverse relay in a nonlinear control 
system. For simplicity in programming, a vertical (or infinite) 
lift curve slope was used locally. This allows the original 
lift curve table look-up data to be used and a simple logic 
statement decreases the interpolated value of lift coefficient 
by an amount AC, if the angle of attack were re than a 
specified amount. Appendix D shows the function subprogram 
Fl and the subroutine TRIM modified from Appendix C, each by 
only one statement. Figure 4 shows the resulting lift curve 
generated by this modification. 

Using the same initial conditions from the previous two 
runs, the same 1/2-degree elevator impulse at three seconds, 
and changing only the lift curve to that in Fig. 4, a limit 


cycle oscillation, or "rocking-chair" stall was generated. 
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The AC, used was 0.050. Figures 5b through 5e show lift 
coefficient, angle of attack, pitch angle and velocity as 
functions of time. The small elevator impulse begins an 
oscillation which grows to a maximum amplitude for the dura- 
tion of the flight time. Figure 5a shows the phase plane 
plot of pitch angle ( © ) versus pitch rate (@). The curve 
winds onto a closed path trajectory, indicative of a limit 
cycle. 

Amplitude of the limit cycle oxcillation could be control- 
led by the magnitude of the AC, on the Lift “eurve? It 
Ac, were reduced by one-half, from 0.050 to 0.025, the 
peak amplitude of the oscillation would be halved. Conversely, 
i aC, were doubled to 0.10, the amplitude of the limit 
cycle would be doubled. The period was not noticeably af- 
fected. Figure 6 shows the time history of the pitch angle 
for varied ac, values. The initial conditions and elevator 
impulse remains the same as for previous runs. 

The magnitude of the elevator impulse had. little effect 
on the limit cycle oscillation. Figure 7 shows the time 
histories of pitch angle for impulses of 0.05 degrees and 2.0 
degrees. The amplitude and period of the limit cycle remain 
the sames The only change was within the first fifteen to 
twenty seconds of flight time when the phugoid mode influenced 
the motion. After the phugoid has sufficiently subsided, the 
limit cycle motion is analagous to that previously encountered. 

To further demonstrate the character of the limit cycle 
and its independence of initial conditions, the aircraft must 


damp to a constant amplitude if the initial oscillation were 
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greater than that for the limit cycle. Figure 8 shows the 
resulting time histories for an aircraft trimmed initially 
at an angle of attack equal to 0.250 radians, well below the 
stall region. An up elevator (back stick) ramp input was 
inserted and held at a new elevator setting corresponding 

to a trim condition close to stall. Since the aircraft seeks 
the trim condition for the elevator angle prescribed, the 
final elevator angle must lie in a region analagous to the 
elevator angle for trim at stall or wherever the instability 
occurs. Because of the rapid change in elevator angle, the 
initial magnitude of the oscillation was approximately three 
times that for the limit cycle. This disturbance also in- 
troduced a large amount of phugoid motion in angle of attack 
and pitch angle as the entire limit cycle was oscillating 
with a long period movement superimposed upon the short 
period behavior. 

A ramp elevator input produced the same Limit cycle. 
Figure 9 shows the time histories for the aircraft initially 
trimmed at an angle of attack equal to 0.400 radians. After 
three seconds, the elevator was deflected at a constant rate 
of one degree per second for 0.30 seconds, and then held 
constant. Again, the final elevator angle is in the region 
of the elevator angle for trim near stall. 

Figure 10 shows the aircraft trimmed at an angle of attack 
well into the stall region ( ®=0.4150 radians). Down ele- 
vator (forward stick) was inserted, by a step input, to an 
elevator angle near that for the stall region. The aircraft 


pitched down and the oscillation grew to;-the,limit,cycle. 
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In summary, the mathematical model used in this analysis 
has created the desired limit cycle situation. The limit 
cycle is shown to be independent of the initial conditions by 
oscillations both decaying and growing to a stable trajectory 
for the state variables. It is self-sustained as shown when 
a small impulse is inserted for 0.05 seconds and then re- 
moved. The phase plane shows a closed path trajectory, 
characteristic of a limit cycle. It must be noted that the 
actual numbers calculated are not of prime importance, rather 
the ability to simulate a known flight characteristic by 
proper aerodynamic nonlinearities was considered valuable. 
The computer program could be used for any aircraft with the 
proper data available, whether from the wind tunnel, cal- 


culated data, or actual flight test. 


IV. CONCLUSIONS 


For a multi-degree of freedom system that is nonlinear 
in its forcing function terms due to position, and has linear 
damping, it is possible to obtain a limit cycle type of 
oscillation in the system's state variables. The governing 
factor in the dynamic stability of the aircraft model was 
the slope of the lift curve. When this slope became -5.0 
or less, the aircraft was dynamically unstable in the region 
OS: Stall. A sharp drop in Lift. coefficient just after 


reaching C produced this effect. This type of stall is 


L max 
characteristic of leading edge flow separation or stall. 
To produce a limit cycle motion or "rocking-chair" type 


of stall, the unstable portion of the lift curve must be 
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bounded by a stable portion on each side. With this type of 
lift curve it is possible to create a limit cycle oscillation 
in the state variables: angle of attack, pitch angle and, to 
some degree, velocity. The frequency and the characteristics 
of this oscillation were analagous to that of short period 
motion, which was predicted by analyzing the roots of the 
equations of motion. 

Since the aircraft is programmed to seek a trim condition 
at the elevator angle specified, the final elevator angle, 
after manipulation, must be close to that for trimmed flight 
in the hector of maximum lift. Various types of elevator 
movement such as impulse, ramp, step, and forward stick 
produced the desired limit cycle. The resulting oscillation 
was independent of the initial conditions imposed, and dem- 
onstrated both a growing and damping oscillation to the 
same amplitude. 

The limit cycle amplitude depended directly on the length 
of the unstable portion of the lift curve. If Ac, were 
doubled, the amplitude would be doubled. This numerical 
experiment provided enough clues to encourage subsequent 
studies for analytic solutions in this area and possibly 
related topics such as post stall-gyrations, lateral and 


directional movement, and spins. 
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